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We consider the long-run growth rate of the average value of a random multiplicative process Xi+i = 
aiXi where the multipliers ai = 1 + pexp(aWi — -a 2 ti ) have Markovian dependence given by the 
exponential of a standard Brownian motion Wi. The average value ( x n ) is given by the grand 
partition function of a one-dimensional lattice gas with two-body linear attractive interactions placed 
in a uniform field. We study the Lyapunov exponent A(p,/3) = limn-^ i log(:En) at fixed /3 = 
| a 2 t n n , and show that it is given by the equation of state of the lattice gas in thermodynamical 
equilibrium. The Lyapunov exponent has discontinuous first derivatives along a curve in the (p, 0) 
plane ending at a critical point (pc, /3c), which is related to a phase transition in the equivalent 
lattice gas. Using the equivalence of the lattice gas with a bosonic system, we obtain the exact 
solution for the equation of state in the thermodynamical limit n —> oo. 


I. INTRODUCTION 

Stochastic recursions of the form aq + 1 = cqaq + 6,; 
with Xi real-valued or matrix-valued quantities and a*, bi 
random variables have been widely considered as mod¬ 
els of dynamics for various processes in physics, ecology, 
computer science, and economics E®. The most stud¬ 
ied case corresponds to i.i.d. random coefficients <q,6,; 
i, but the case with state dependence (for example of 
Markovian type) has been also considered Such pro¬ 
cesses can produce a wide variety of distributional prop¬ 
erties for x„ and under certain conditions on the coeffi¬ 
cients cq, bi they can generate heavy-tailed distributions 
mm. Alternative mechanisms for generating heavy 
tailed distributions with specific application to financial 
time series are discussed in 0,1- 

We consider in this paper the discrete time random 
multiplicative process defined by 

x i+ i = ctiXi , cq = 1 + pe aWi ~i a2ii (1) 

with Wi = W(li) the values of a standard Brownian 
motion starting at W(0) = 0 and sampled at uniformly 
spaced times A = it. The parameter p is positive and 
bounded as 0 < p < 1. This is a random multiplicative 
process with non-stationary multipliers cq, which have 
Markovian dependence introduced through the depen¬ 
dence on the Brownian motion W(t). 

The model o can be used to describe the discrete 
time dynamics of a quantity which changes in each pe¬ 
riod by an amount which is proportional to its value 
at the beginning of the period. The random multiplier 
is positive and follows a geometric Brownian motion in 
discrete time. The simplest such variable is a bank ac¬ 
count which accrues interest by simple compounding over 
each period (tj,tj + i), assuming that the interest rate 
Lj,i +1 = Lq exp(aWi — ^cr 2 U) follows a geometric Brow¬ 
nian motion in discrete time. This corresponds to the 
Black, Derman, Toy model of stochastic interest rates in 
mathematical finance [16]. More generally, the model 0 
could be used to describe multiplicative processes with 


Markovian dependence, such as for example of the type 
considered in 11] in the context of population dynamics. 

The random multipliers cq in Eq. o have the same 
average value (cq) = 1 + p. If the multipliers ai were inde¬ 
pendent random variables, one would expect the average 
of the process 0 to have a simple geometric growth as 
(x n ) = Xo(l + p) n - However, the successive multipliers cq 
are correlated through their dependence on the common 
Brownian motion W,. Their covariance is 

co v(ai,a,j) = (didj) - (a,i){aj) (2) 

= p 2 (e a2min{ti ’ tj ) - 1) 

The average value (x n ) of the process © can be com¬ 
puted exactly, and displays a surprising explosive behav¬ 
ior for sufficiently large volatility a or sufficiently large 
time step n [THI ] . This is very different from the naive 
expectation of a geometric growth obtained by assum¬ 
ing statistically independent multipliers, and is an effect 
due to the correlation between successive multipliers cq. 
The same explosive behavior is noted also for the higher 
positive integer moments (( x n ) p ). This implies that the 
probability distribution of the variable x n develops heavy 
tails under certain conditions. 

The process 0 is somewhat similar to the Kesten 
process 0], which is defined by the stochastic recursion 
x ^|_i = diXi + bi with di,bi i.i.d. random variables. For 
this case it has been shown a that under certain con¬ 
ditions on the distributions of cq ,bi, the distribution of x n 
approaches a stationary form and has heavy tails (of reg¬ 
ular variation). This class of models has been extended 
to accomodate also Markovian dependence for the cq,&i 
factors, see e.g. 0 m. 

In this paper we consider the long-run asymptotics of 
the average value (x n ) of the process 0 in the limit 
n —> oo while keeping j3 = \a 2 t n n constant. This is de¬ 
scribed in terms of a Lyapunov exponent A (p, /?), defined 
below in Eq. ©. We study the properties of the Lya¬ 
punov exponent and derive an exact solution for \(p, (3) 
in terms of the equation of state of an equivalent one- 
dimensional lattice gas. The Lyapunov exponent is a 
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continuous function of its arguments and has a discon¬ 
tinuous derivative along a curve in the (p, (3) plane ending 
at a critical point. The qualitative features of the exact 
solution are well reproduced in terms of an approxima¬ 
tive mean-field theory with van der Waals equation of 
state. 


II. LYAPUNOV EXPONENT AND SCALING 

We study in this paper the long-run growth rate of 
the average value ( x n ) of the process © as the number 
of time steps becomes very large n —> oo. This can be 
described in terms of a Lyapunov exponent A, defined as 

A = lim - log(i n ). (3) 

n —>-00 fi 

A similar quantity and its properties were studied in 
the context of random multiplicative process with i.i.d. 
matrix-valued multipliers ctj in (ToL fl3l. [T3| and for pro¬ 
cesses with Markovian dependence in |ll|. 

For the study of the n —> 00 limit it will prove helpful 
to make use of the equivalence between the process © 
and a lattice gas noted in |15] . We recall briefly the main 
points of this analogy. Consider a one-dimensional lattice 
gas with n — 1 sites, and denote the occupation number 
of the site i as m = {0,1}. The lattice gas particles 
interact by the Hamiltonian H = Tiirij£ij where 

the interaction energy of two particles at sites i,j is 

= (4) 

3 { +00 , * = j w 

Writing — min(i,j) = \\i — j\ — i+j ) the Hamiltonian 
can be put into the form 

H = ^ n i n 3Vij + J! ( 5 ) 

i<j i 

with translation invariant 2-body interaction ipij = | i— 

j | and single-site energies pi = —The lattice gas 
particles attract each other with linear two-body inter¬ 
actions and are placed in a uniform strength field which 
drives them towards the right side of the lattice. 

The precise relation of the random multiplicative pro¬ 
cess © to the one dimensional lattice gas with interac¬ 
tion © is given by the equality 

(x n ) = x 0 Z(/3,p ), (6) 

between the average ( x n ) and the grand partition func¬ 
tion Z(/3,p) of the lattice gas with fugacity p and tem¬ 
perature T = 1/(3 with 

P = ^ 2 t n n. ( 7 ) 

The grand partition function of the lattice gas is 

n— 1 

m P ) = T / p N z N ((3), (8) 

AT =0 


with Zn(/3) the canonical partition function, given by a 
sum over all configurations with N occupied sites 

Zn(J 3)= e ~ PH - ( 9 ) 

{«i}>XT n i =N 

The relation © can be used to prove the following 
result: The Lyapunov exponent, defined as the limit 

A (p,T)= lim - log (x n ) = l-p(p 7 T) (10) 

n—^oo n 1 

exists and depends only on p and T = 1/(3 with (3 given 
by ©• Furthermore, A(p, T) is expressed as shown in 
terms of the pressure of the lattice gas p(p, T). 

This result follows from the existence of the thermo¬ 
dynamical limit for the lattice gas with interaction ©. 
A sufficient condition for the existence of this limit for a 
lattice gas with translation invariant 2-body interaction 
Pij is the finiteness of the sum 

XllV’bl <°°, ( n ) 

for any site i. One can easily check that this condition is 
indeed satisfied by the interaction ©. The factor 1/n 2 
was introduced in Eq. © motivated by this condition. 
Under this condition, the limit 

lim - log Z(/3, p) = (3p(p, T) (12) 

n—>oo fi 

exists and is a finite function of temperature T and fu¬ 
gacity p 0- This function is related to the pressure 
of the lattice gas p(p,T), which is given in terms of the 
thermodynamical potential as 

ft =-T log Z = -np(p,T). (13) 

The existence of the limit in m together with © im¬ 
plies the result (fTUl) . It is clear also from this derivation 
that the limit in (flUl) does not depend on the initial value 
xq of the random multiplicative process. 

The result (fTUl) shows that there is a close relationship 
between the long run asymptotics of the random multi¬ 
plicative process © and the equilibrium thermodynam¬ 
ical properties of the lattice gas with interaction © in 
the thermodynamical limit n —> 00 . It also implies that 
the Lyapunov exponent A has a scaling property in the 
large n limit, as it depends on cr, r, n only through the 
combination (3 = 1/T defined in ©. 

Numerical simulations show that the average (x n ) has 
an explosive behavior for sufficiently large cr or n [15| |. We 
will show here that this phenomenon can be related to 
non-analyticity of the Lyapunov exponent A(p, T) in its 
arguments. In the next section we summarize the results 
of the numerical simulations, and in Sec. lIVI we present an 
analytical approximation based on a lower bound for the 
partition function of the lattice gas which corresponds 
to a van der Waals system, which reproduces the qual¬ 
itative features of the numerical simulation. In Sec. m 
we present the exact solution for the Lyapunov exponent 
A (p, T) following from the equation of state of the lattice 
gas with interaction © in the thermodynamical limit. 
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III. NUMERICAL RESULTS 

The integer positive moments (( x n ) p ) of the random 
variable x n defined by the process © can be computed 
exactly, but numerically, using a recursion relation de¬ 
scribed in Appendix A of [15;]. We use this method to 
compute the average value (x n ) and the finite n approx¬ 
imation to the Lyapunov exponent X n = — log(:r n ). 

We show in Figure [T] (upper panel) typical numerical 
results for A„ as function of temperature T = 1//3 defined 
as in © for several values of p between 0.005 and 0.125. 
The simulation has n = 200 time steps of size r = 1 and 
initial value xo = 1. The results of the simulation show 
that A„ is always positive and is a decreasing function of 
T, approaching a small but finite value as T —> oo. The 
functional dependence is qualitatively different, depend¬ 
ing on whether p is below or above some critical value 
pc- This suggests the following picture: 

i) p > pc- The Lyapunov exponent A(p, T) is a smooth 
decreasing function of temperature. 

ii) p < pc- The T dependence of A(p, T) has a kink at a 
certain transition temperature T tl -(p), and its derivative 
with respect to T is discontinuous at this point. As T 
decreases below this value, A(p, T) increases very rapidly 
and explodes to infinity as T — > 0. 

The approach to the thermodynamical limit n — > oo 
is illustrated in Fig. |T] (lower panel) which shows log A„ 
vs the temperature T for p = 0.025 as n is increased: 40 
(black), 100 (blue), 200 (red). These plots confirm the 
scaling property of the Lyapunov exponent, more pre¬ 
cisely that A„ —> A(p, T) approaches a function which 
depends only on (p, T) in the n —> oo limit. The numer¬ 
ical results show that with n = 200 the scaling property 
holds very well. 

This picture can be used to understand the behavior 
of the expectation value (x n ) observed in numerical sim¬ 
ulations as functions of the volatility a, maturity t n and 
time step r [r|. As the product (corresponding 

to the scaling variable /3 defined in Eq. © which has the 
meaning of the inverse temperature) increases, the tem¬ 
perature of the equivalent lattice gas decreases. As long 
as the temperature T is above the transition tempera¬ 
ture T't r (p), the lattice gas is in the gaseous phase and 
the pressure increases but remains small. The equivalent 
statement for the process © is that the average value 
(x n ) increases slowly. As soon as the temperature drops 
below the transition temperature, the lattice gas con¬ 
denses into the liquid phase and the pressure increases 
very fast with the temperature, which translates into an 
explosive growth of the expectation (x n ). The approach 
to the n —» oo limit of the random multiplicative process 
© is discussed in more detail below in Sec. ED 


IV. VAN DER WAALS APPROXIMATION 

We present in this Section analytical approximations 
for the thermodynamical properties of the lattice gas 




T 


FIG. 1: Above: The Lyapunov coefficient X n = A log^™) as 
function of temperature T for several values of the p param¬ 
eter p = 0.005,0.0125,0.025,0.05,0.125 (from left to right). 
The simulation has n = 200 time steps with r = 1 and xo = 1. 
Below: the approach to the thermodynamical limit. Plots of 
log An vs T at p = 0.025 for n = 40 (black), n = 100 (blue), 
n = 200 (red). 


with long-range interactions 0 - These approximations 
are based on lower and upper bounds on the partition 
function of the lattice gas with the interaction 0 

Z(fi,p)>Z{p,p)>Z{p,p). (14) 


These bounds are partition functions for lattice gases 
with uniform infinite range interactions. The lower 
bound has interaction jl5| 


/ 3n ^ 3 
\ +oo ,i = j 


(15) 


The lower bound on the energy of a state with N particles 
En > = gi T V(V — 1)(2V + 2 — 3n) > gives an 

upper bound on the canonical partition function Zn{/3) 
in terms of the partition function of a lattice gas with 
uniform interactions 

~ 3 F 2 — 2 N — 2) ,i^j 

Too ,i = j ' 

A weaker upper bound which has the advantage that it 
does not depend on the number of particles (similar to 
the lower bound (ll5l) ) is ^ j. 
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Each of these simpler systems is equivalent to a mean- 
field theory, as each particle feels the effect of the other 
particles as a constant interaction energy. Their ther¬ 
modynamical properties are given by the van der Waals 
theory 01 H- The lower bound is saturated in the very 
large temperature limit T —> oo [HI, and is numerically 
more accurate for all temperatures. For this reason we 
will restrict ourselves to the predictions following from 
the lower bound. 

The thermodynamical quantities of the approximative 
van der Waals system can be computed in closed form. 
The free energy corresponding to the lower bound (fl5l) 
is F(n,N,T) = nf(d , T) with f(d,T) the free energy 
density, given by (T5|, L 3 J]| 

f(d,T) = CE[T(dlogd+ (1 - d) log(l - d)) - id 2 ](17) 

with d = N/n the lattice gas density and CE[f(d,T)\ 
denotes the convex envelope of f(d , T) with respect to d. 
The equation of state has van der Waals form 

P = —T log(l - d) - ^d 2 , (18) 

supplemented with the Maxwell construction, which re¬ 
places (fT8l) with a constant pressure po(T) in the region 
d g < d < de with 

d g {T) = \{ 1-A), de(T) = ^(1 + A) (19) 
where A is the positive solution of the equation A = 



The fugacity p is 


P = 


1-d 


( 20 ) 


The inversion of this relation in order to find d for given 
(p, T) requires some care. For T > Tq = 1/6 this equa¬ 
tion has a unique solution for d. For T < 1/6 it has 
three solutions. Denoting the smallest and largest solu¬ 
tions with di(p,T) and d 2 (p,T), respectively, the lattice 
gas density is given by 


d(p,T) 


di{p,T) ,p>e 3 t 
di(p,T) , p < e“3T 


( 21 ) 


At p = e 3 t the gas and liquid phases can co-exist, and 
their densities are 

d 2 [e~^ ,T) = d t (T) , d x {e~™ ,T) = d g (T). (22) 


These results can be used to understand the qualitative 
behavior of the curves for the Lyapunov exponent A(p, T ) 
in Fig. [I] In the van der Waals approximation this is 
given by 


AvdW 


-log(l-d)-^, d<£ {d g (T),d e {T)) (2 „, 

Apo(T), (d g (T),d e {T)) 1 J 


with d = d(p,T) given by Eq. (1211) . and po(T) corre¬ 
sponds to the flat portion of the isothermal curves and 
is determined such that A v dw is a continuous function of 
d. The shape of the curves A v dW vs d is shown in Fig. [2] 
(lower plot) and are essentially the well-known van der 
Waals isothermal curves. 

The function \ v dw(Pi T ) at fixed p has a discontinuous 
derivative with respect to T at the transition temperature 


m(vdW) 

J tr 


(P) = 


1 

31ogp ’ 


(24) 


provided that the fugacity is below the critical value p < 
pc = e~ 2 . Alternatively, at fixed T the transition occurs 
at the point p = exp(— ^) where the solution for d in 
m changes branches and switches between d\ and d 2 . 
At this point the pressure is a continuous function of T 
but its derivative has a jump. The maximal value of p 
for which the transition occurs corresponds to the critical 
point of the van der Waals system which has parameters 

Tc = g , d c = -■ pc = e -2 ■ (25) 

These results are illustrated in Fig. [2] The upper plot 
shows log A (p,T) as function of the temperature T at 
fixed fugacity p. The solid lines correspond to the nu¬ 
merical simulation for the exact lattice gas with interac¬ 
tion Q while the dashed curves correspond to the van 
der Waals approximation (1231) . Both these curves show 
a sharp transition at an intermediate transition temper¬ 
ature T tI (p) which depends on the fugacity p. 

In order to understand better the discontinuous be¬ 
havior of the pressure p(p, T ) with respect to T at fixed 
fugacity p we show in the lower panel of Fig. [2] plots of 
A v dw vs the density of the lattice gas d for several val¬ 
ues of the temperature T = 0.08 — 0.15 in steps of 0.01 
(from bottom to top). The top-most curve corresponds to 
T = Tc = 1/6. The colored dots are on the curve of con¬ 
stant fugacity p = 0.05 and correspond to increasing tem¬ 
perature. The points to the right of the horizontal por¬ 
tion correspond to T = 0.105,0.11,0.1105,0.1108,0.111 
and the points to the left correspond to T = 0.12,0.3 
(from top to bottom). The transition temperature for 
this fugacity is T^ dW \ 0.05) = —l/(31og0.05) = 0.1113. 
As T increases above the transition temperature, the den¬ 
sity drops suddenly from the right side (liquid phase) to 
the left side of the flat portion of the isothermal curve 
(gas phase), and the lattice gas evaporates. As explained 
above, this jump is responsible for the discontinuity in 
the slope of the curve p(p, T ) vs T at fixed p, occuring at 
the transition temperature Ttr(p)- 

Using the van der Waals picture we can obtain analyti¬ 
cal expressions for the Lyapunov exponent A v dw (T, p) in 
the small- and large-temperature limits, which describe 
the right and left tails of the curves in Fig. [2] (upper plot). 

In the large-temperature limit T —> oo the lattice gas is 
in the gas phase. From Eq. (l23l) the Lyapunov exponent 
becomes approximatively equal to A v dw — — log(l — d) 
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FIG. 2: Above: Plots of log A vs the temperature T at fixed 
fugacity p. Solid curves: numerical solution for the lattice 
gas with interaction 0 with n = 200 sites, dashed curves: 
the van der Waals approximation A v aw corresponding to the 
system USD- The colors correspond to the same values of the 
fugacity p as in Fig. [T] Below: The Lyapunov exponent A v dw 
vs the lattice gas density d at fixed temperature. The curves 
shown correspond to T = 0.08 — 0.15 in steps of 0.01 and 
Tc = 1/6 (from bottom to top). The colored dots are on 
the p = 0.05 curve, see explanation in text. The dashed line 
shows the ideal gas approximation A = — log(l — d). 

in this limit. Inversion of Eq. (BUI) gives that the density 
approaches a constant value d(p,T ) = + 0(T _1 ). 

Substituting into BUI) gives the large temperature limit 
of the the van der Waals approximation of the Lyapunov 
exponent 

lim A v dw(p, T) = log(l + p). (26) 

T—H-oo 

The independence on temperature of this limiting value is 
due to the fact that p/T approaches a universal function 
in the small density limit, which is just ideal gas behavior. 
This is seen as the overlap of all the isothermal curves in 
the low density region d «C 1 in the lower panel of Fig. [2] 
All these curves approach the ideal gas equation of state 
p/T = d which is shown as the dashed line. 

The result (BUI) agrees with the a = 0 limiting behavior 
of the random multiplicative process 0 which is solved 
trivially as x n = Xo(l + p) n . This gives Al^o = log(l + 
P)- 

In the small-temperature limit T —> 0 the lattice gas 


is in the liquid phase. The density d is close to 1, and is 
given by the approximative solution of BUI) 

hm d(p, T) = 1 - 1 2 ■ (27) 

t->o l + pe3 p 

Susbtituting into the expression BUI) we obtain the small- 
temperature asymptotics of the Lyapunov exponent 

lirn A vdw (/J, T) = log(l + pel 8 ) - ~ logp+ -^(.28) 

T—>o d 61 

This gives an explosive behavior of A vd w {p, T) as T —> 0 
which reproduces the results of the numerical simulation 
shown in Fig. [2] (upper plot). 

The relations m and BUI) give the small- and large- 
temperature asymptotic behavior of the van der Waals 
approximation for the Lyapunov exponent of the random 
multiplicative process. As seen from Fig. [2] (upper plot), 
they reproduce reasonably well the main qualitative fea¬ 
tures of the curves obtained from the exact numerical 
simulation of the model. 


V. EXACT SOLUTION IN THE 
THERMODYNAMICAL LIMIT 

The thermodynamical properties of the lattice gas with 
interaction (0) can be found exactly in the thermodynam¬ 
ical limit n, N — > oo at fixed particle density d = N/n. 
We start by stating the solution. 

Proposition 1. Define the intensive quantity n(d 1 T) 
as the solution of the equation 

i _ ! = f 1 _ d JL _ (29 ) 

d J 0 e^ d2 y^ 2 ~y)+' K \ — l ’ ^ 

with d = N/n the density of the lattice gas. 

The free energy density is given by the convex envelope 
with respect to d of the function 

f(d,T) = ^d 2 (2d-3) + n(d-l) (30) 

+ Td [ dy\og(l ^ e-^ d2y ( 2 ~ y)+7r] ), 

Jo 

and the pressure is 

p(d,T) = -f(d,T) + dd d f(d,T) (31) 

= -d 2 (M-3) + n + 2d 3 ^ e/3[d 2 ^ 2 - y )-H _ ' 

The Gibbs free energy density g(d,T) defined as G = 
E — TS + nn = ng(d , T) = n(f(d, T) + n) is given by 

g(d, T) = \d 2 (2d-3)+Trd (32) 

+ Td f dy \og(l - e - p[d2y{2 - v)+ ^). 

Jo 
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The chemical potential y = ddf{d,T) is given by the 
change of the free energy density at fixed volume n and 
temperature T when adding one particle. The lattice gas 
fugacity is 

P(d , T) = = exp (±d d /(d, T)) . (33) 

We will derive these results using the isobaric- 
isothermal ensemble, see for example [20]. The proof is 
based on the equivalence of the lattice gas with n— 1 sites 
and TV particles with a boson system of n — TV — 1 bosons 
which can be placed onto TV + 1 known energy levels 
ujk,k = 0,1,--- ,7V. The isobaric-isothermal ensemble 
for the lattice gas will be shown to be equivalent with 
the grand canonical ensemble for the equivalent bosonic 
system. 

We start by recalling the energy spectrum of the lat¬ 
tice gas with interaction a. The energy of the lattice 
gas with interaction a and TV particles is given by (see 
Proposition 1 of 0) 


N 

E n = + J2 VkUk (34) 

fc=o 

where yk = 0,1, 2, • • • are the occupation numbers of the 
energy levels cok which are given by 

uj k = \k{N- bfc + 1)), & = 0,1, ■ • • , TV. (35) 

2 

The ground state energy is 

E^ = ^ N { N ~l)(2 N + 2-3n). (36) 

The occupation number yk of the k —th energy level has 
a geometrical interpretation as the number of the empty 
sites between the (k — 1)—th and k —th particles on the 
lattice (ordered in increasing order i\ < i^ < ■ ■ ■ < In)- 
(For k = 0 the variable yo has the meaning of the number 
of empty lattice sites to the left of the leftmost particle.) 
They satisfy the sum rule J/fc = n — N—1 which has 

a simple geometrical interpretation as the total number 
of the empty sites of the lattice gas. It is clear that the 
states and the energy of this system are exactly equiva¬ 
lent to those of a system of n — TV — 1 bosons which can 
be placed onto the TV + 1 energy levels uik- 

In the isothermal-isobaric ensemble one fixes the num¬ 
ber of particles TV and the temperature T, while the lat¬ 
tice size n (volume) is allowed to be variable, subject 
only to the constraint n > TV. Furthermore, an intensive 
quantity tt is introduced, which is fixed and is usually 
identified with the pressure of the system 1 . One defines 


1 In our case, due to the volume dependence of the interaction Q, 

this quantity will be seen to be different from the pressure, so we 
denote it with a different symbol. 


the isothermal-isobaric partition function as 

oo 

5(t r,T) = Y, Z N (m)e- 0 ™ . (37) 

m=N 

Here Z^fn) is the usual canonical partition function of 
the lattice gas with n sites and TV particles. The Gibbs 
free energy G and the free energy F are determined as 

G = F + Tr(d, T)n = —T log S(7T, T) (38) 

where the intensive quantity 7r(d, T) is determined from 
the equation for the average lattice size 


n= -a w 5(7r,T). (39) 

In order to be able to perform the sum over m in (1571) 
in closed form we replace m — > n in the denominator of 
a, where n is a fixed value which will be set equal to the 
actual lattice size n. As a result the same replacement is 
made in the denominators of (1551) and (1551) . We will iden¬ 
tify in the intermediate steps N/n = d with the density of 
the lattice gas. The replacement of the summation index 
m with the average value of the lattice size n is justified 
in the thermodynamical limit, when the fluctuations of 
the volume m around its mean become vanishingly small. 
However, this replacement has the effect that the inten¬ 
sive quantity 7r is not exactly equal to the pressure of the 
lattice gas. 

Using the representation (1551) we write the canonical 
partition function as 

Zn{iti) = e _/3B «‘ > Zjy(m) (40) 

= exp ^ - 7^2 7V(7V - 1)(27V + 2 - 3m)^ Z N (m ), 

where we introduced Z^(rn) the partition function of a 

system with energy E N = J2k=o V^k and J2k=o V* = 
m — TV — 1 . 

The sum over m in (EZD can be evaluated in closed 
form using the combinatorial identity 


Y, Z N (m)e~^ 

m=N 


g-flNn 

n£L 0 (i - e-h^+'K)) ‘ 


(41) 


This has a clear resemblance to the grand partition func¬ 
tion of a system of non-interacting bosons which can be 
placed onto TV + 1 energy levels ojk, see e.g. [21] . The 
intensive quantity n is the analog of the (minus) chemical 
potential for the boson system. 

The isothermal-isobaric partition function of the lat¬ 
tice gas is 


S(7t,T) = e _ 3( 3d2 ( Ar+1 ) . 


0 —pN (-K—d 2 ) 


(42) 


n£L 0 (i - e -0G k +K-d^ 

The Gibbs free energy is obtained as 

G(tt,T) = -TlogH(7r,T) (43) 





7 



0.2 0.4 0.6 0.8 1 

d 



0.2 0.4 0.6 0.8 1 

d 



FIG. 3: The exact isothermal curves p(d, T ) (black) of the lat¬ 
tice gas © obtained from Eq. GD- The horizontal black line 
in the upper plot was obtained by the Maxwell construction 
as described in Sec. IV Bl The numerical solution on a lattice 
gas with n = 200 sites is shown as the red curve. The van 
der Waals approximation of the isothermal curves p v dw(d, T) 
given in Eq. m is shown as the dashed blue curves. The 
three plots correspond to T = 0.1, 0.2, 0.3. 


= \d 2 (N+l) + N(n-d 2 ) 

O 

N 

+T^log(l - e -/3(^+^-rf 2 )). 

k=0 

In the N —»• oo limit the sum can be replaced with an 
integral as 

N 

J2 l°g(l - e-^ k+n - d2 )) ~ (44) 

k=0 

N dy log (l - e - 0 ld 2 y{ 2 -y)+n-d^ _ 


Collecting together all the terms we find the result for 
the free energy density 

f{d,T) = —-n =^d 3 +7r(d-l) (45) 

n 3 

+ Td J o dy\og(l-e-^ d2y ^- y)+7T - d2) y 

The convex envelope of this function with respect to d is 
understood. This reproduces the relation Eq. (1501) upon 
replacing 7r — d 2 —> tt. Since 7r is not an observable quan¬ 
tity we can redefine it by a shift tt — d 2 —> n, in order to 
simplify the form of the result. Alternatively, the result 
Eq. m can be directly obtained by replacing m —> n 
also in the numerator of in Eq. (l40l) . 

The pressure is obtained as 

p(d,T) = -f(d,T) + dd d f(d,T). (46) 

Substituting here f(d,T) from Eq. fi5l) and replacing 
again n — d 2 —> tt gives the result (1511) . 

The intensive quantity tt is determined from the con¬ 
dition that the average lattice size is equal to its actual 
value (m) = n. This is expressed as 

n = ~ Td - lo § T ) = N + E e P(u, k+ Ld 2 ) _ ! • ( 4? ) 

k =0 

This has again a clear resemblance to the Bose-Einstein 
distribution, as it expresses the total number of bosons 
n — N as a sum over the Bose-Einstein distribution for 
all energy levels Wfc 


N = e 0(u k +lT~d 2 ) _ l ■ 

k =0 

Replacing the sum with an integral in the N —> oo limit 
gives 


l 

n — N = N / dy a ,,, , —r--- (49) 

J y e P(d 2 y{2-y)+Tv-d 2 ) _ i \ > 

which reproduces the equation (1311) with the substitution 
7 r — d? —y tt. This completes the proof of Proposition 1. 


A. Properties of the exact solution 


We prove here a few properties of the exact solution 
presented in Proposition 1. 

The equation 


i _ i = r __ 

d J Q e P[d 2 y{2-y)+Tt] _ i ■ 


(50) 


has a unique solution 7r(d, T) > 0 for any 0 < d < 1 and 
T > 0. This implies that there is no analog of the Bose- 
Einstein condensation in the equivalent bosonic system. 
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FIG. 4: The exact fugacity p(d, T) of the lattice gas (black 
curve) from Eq. (1331) . The horizontal black line in the upper 
plot was obtained by the Maxwell construction described in 
text. The numerical solution for p(d,T) for a lattice gas with 
n — 200 sites is shown as the red curve. The van der Waals 
approximation p(d,T) given in Eq. (1201) is shown as the dashed 
blue curves. The plots correspond to T = 0.1, 0.2. 


Furthermore, the solution 7 r(d, T ) is bounded from below 
as 


7 r(d, T) > n(d, T ) = — T log(l — d) — -cP . (51) 

o 

To prove this result, consider the integral 

= j Q e f)[d?y{2-y)+n] _ 1 ■ ( 52 ) 

This is a monotonously decreasing function of 7 r for 7 r > 
0 , which approaches zero as 7 r —> oo and diverges to +oo 
as 7 r —> 0. This implies that the equation 1/d — 1 = /( n) 
has a unique positive solution for 7 r for any value of d € 
( 0 , 1 ). 

The function I (tv) is bounded from below as 


/( 7 r) > I(tt) = 


1 


e P[%d 2 +n] _ 


(53) 


This follows from the Jensen inequality, as 
( e /?[d y(2-y)+K] _ i)-i j s a convex function of y. 
As a result we have 


dy 


> 


e 0[d 2 y(2-y)+n] _ l — e /3[§d 2 -|-7r] _ ' 


(54) 


Using this into (1311) gives the inequality (l5U) . 

Numerical simulation shows that the inequality (1511) 
approaches saturation in the large temperature T —> oo 
limit. 

Define the integral 

J(x) = J dy log (l - e -Pd 2 y(i-y)+^ . ( 55 ) 

Numerical simulation shows that in the large tempera¬ 
ture limit the integral (1551) approaches the limiting value 

lim J(f3n) = logd. (56) 

T —>oo 

Substituting 7 r —► 7 r and J logd, the free energy 
density becomes 

lim f(d,T) = TA(d)-U 2 (57) 

T—> oo d 

with A{d) = dlogd + (1 — d) log(l — d). This coincides 
with the van der Waals result corresponding to the lower 
bound on the canonical partition function (1151) . This 
bound implies that the free energy density is bounded 
from above as 

f(d,t) < TA(d) — . (58) 

The bound is saturated in the large temperature limit 
and gives the van der Waals equation of state of the lat¬ 
tice gas 0 

P vdw(rf,T) = -Tlog(l-d)-id 2 . (59) 

It is possible that the inequality (l58l) can be proved also 
analytically starting from the definition of ir(d, T) given 
by (13U) but the author was unable to do so. However, the 
numerical solution confirms that this inequality is indeed 
satisfied in all cases considered. 


B. Numerical results for the exact solution 

We present in this section the numerical results for the 
thermodynamical quantities of the lattice gas following 
from the exact solution in the thermodynamical limit. 

Fig-EHshows numerical results for the isothermal curves 
p = p(d, T) for several values of the temperature T ob¬ 
tained by solving the equation m■ The curves p(d , T) 
have a qualitative resemblance to the van der Waals 
isothermal curves. For sufficiently small temperature, be¬ 
low a critical temperature T < Tq , the pressure has two 
extremal points as function of d (van der Waals loops), 
while for T > Tc the pressure is a monotonously increas¬ 
ing function of d. The critical isothermal curve T = Tc 
has an inflexion point where the second derivative van¬ 
ishes d%p(d,Tc) = 0. 

The isothermal curves for T < Tc must be supple¬ 
mented by the Maxwell construction. This determines 
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FIG. 5: The phase co-existence curve T tr (p) following from 
the exact solution of the model (black curve). This ends 
at the critical point C with coordinates given in (1611) . The 
blue dashed curve shows the van der Waals approximation 
T t ( r vdw) (p) =-l/(31ogp). 


the flat portion of the pressure curve Po(T) and the gas, 
liquid densities d g (T),de(T) from the equal-area condi¬ 
tion 


r de{T) dx 

/ ™(p(x,T)-p 0 (T))=0. (60) 

Jd g (T) x 2 

The results for d g (T),de(T) can be used to find the fu- 
gacity po(T) at which the two phases are in equilibrium 
at temperature T. The result is shown in Fig. [2 where 
we compare the exact result for p(d, T) (black curves) 
against a numerical solution for a lattice gas with n = 200 
curves (red curves) and the van der Waals approximation 
(dashed blue curves). 

The curve po{T) defines the phase co-existence curve. 
This is shown in Fig. [5] as the black solid curve, together 
with its van der Waals approximation T^ dW \T) = 
— 3 (blue dashed curve). It starts at origin and ends 
at the critical point (pc,Tc). The exact solution of the 
lattice gas given by the Proposition 1 gives the following 
approximative values for the critical parameters 


with c a constant independent of n and A (p, T) depends 
on the two parameters of the equivalent lattice gas sys¬ 
tem: p and T = 1//3 with f} defined in J7J). 


A. Long-run limit n —> oo at fixed time step r 


Consider the behavior of (x n ) at fixed model parame¬ 
ters <r, p, t as n —> oo. This corresponds to the scaling 


p = fixed 
2 

T = 


(63) 

(64) 


As n is increased, the temperature of the equivalent lat¬ 
tice gas decreases. The state of the system is described by 
a point which describes a linear trajectory in the (p, T) 
plane, given by a vertical line ending at T = 0, see Fig. [ 6 ] 
We distinguish two possibilities for the behavior of 
the average (x n ) as n —> oo. If p < pc (the orange 
line marked (I)), the trajectory intersects the phase co¬ 
existence curve, and the Lyapunov exponent A will have 
a discontinuous first derivative at the transition temper¬ 
ature T tr (p) as shown in Fig. [T] At this point the average 
( x n ) has a fast explosive growth with n. If p > pc then 
the Lyapunov exponent has a smooth dependence on T, 
and the average (x n ) increases smoothly with n. 


B. Fixed maturity t and volatility a 


We consider here the asymptotic behavior of (x n ) at 
fixed maturity t n = t and volatility a as the number of 
time steps is increased n —> oo. The size of the time 
step decreases and approaches zero as r = t/n —> 0. We 
distinguish two possible cases: 

i) fixed p. The equivalent lattice temperature scales 
with n as 


T = 


2 

cr 2 tn 


£l 

n 


(65) 


Tc = 0.195, dc = 0.36 , pc = 0.122. (61) 

These are close to the critical parameters of the van der 
Waals approximative system Tc = \jdc = \,Pc = e -2 . 


VI. LONG-RUN AND CONTINUOUS TIME 
LIMITS 

We consider in this section the implications of the re¬ 
sults for the Lyapunov exponent on the asymptotics of 
the average value (x n ) of the random multiplicative pro¬ 
cess O in the long-run and continuous time limits. The 
long-run asymptotics of the average has the general form 

lim 

n—> oo 


and thus the system cools down as the number of time 
steps n increases. The state of the system is described by 
a point which describes a linear trajectory in the (p, T) 
plane, given by a vertical line ending at T = 0, see Fig. [6] 
The behavior of the system is very similar to the case 
discussed above in Sec. IVI Al 

If P < PC (the orange line marked (I)) then the Lya¬ 
punov exponent has a kink at the transition temperature 
T tr (p) as shown in Fig. |T] If p > pc then the Lyapunov 
exponent has a smooth increase with n, without any dis¬ 
continuous behavior. 

ii) scaling p = rr with fixed r > 0. This corresponds 
to scaling both the temperature and fugacity to zero in 
a correlated manner such that their ratio stays constant 

/3 = Cin, p = C 2 — , with (3p = \a 2 t 2 r . ( 66 ) 

n 2 


Mp,T)n 
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FIG. 6: Paths in the plane (p, T) corresponding to the dif¬ 
ferent approaches to the n —> oo limit discussed in Sec El 
The black solid curve shows the phase co-existence curve of 
Fig. 0 and ends at the critical point C. i) orange curve (I): 
fixed time step t, increasing time, ii) dashed blue and red 
curves (II), (III): fixed maturity t n and volatility, decreasing 
time step with scaling p = rr. 


The state of the system is described by a trajectory in 
the (/?, T) plane which is a straight line joining the origin 
with the point (rt, 2/{cr 2 t)). The origin corresponds to 
the continuous time limit n —> oo. 

The qualitative behavior of the Lyapunov exponent 
as n is taken to be very large is different depending on 
whether this line intersects the phase co-existence curve 
or not. These two cases are illustrated in Fig.[G]by the two 
dashed curves marked as (II) and (III). For the case (II) 
the Lyapunov exponent has a discontinuous derivative 
at the value of n which corresponds to the intersection 
with the phase co-existence curve, while in case (III) the 
Lyapunov exponent has a smooth dependence and (x n ) 
increases smoothly as n is increased. 


C. Continuous time limit 

The n —> oo limiting form of the probability distribu¬ 
tion function of x n for case ii) of Sec. lVIBl can be found 
in closed form. With the scaling p = rr, the random 
multiplicative process m approaches in the continuous 
time limit r —> 0 the diffusion defined by the stochastic 
differential equation 

dx(t) = re aWW ~^ 2t x{t)dt . (67) 

Conversely, the random multiplicative process m can be 
regarded as an Euler discretization in time of the contin¬ 
uous time process (EH- 

The solution of (1671) is given by the exponential of the 
time integral of the geometric Brownian motion 


The distributional properties of this quantity are well 
studied in mathematical finance, see (H, [H, HU HU and 
references cited. We summarize here the main results. 

Denote the probability distribution function of the 
time integral of the geometric Brownian motion 

/■t 

A(t) = / dse "W'M-** 2 . ( 69 ) 

J o 

as <F(z,f) = P(A(f) G (z, z + dz)). This distribution 
approaches a simple form in the very small and very large 
time limits respectively [23j . In the small time limit the 
distribution is log-normal with mean t. and variance ^ er 2 f, 
while in the infinite time limit it approaches a stationary 
distribution given by the inverse Gamma distribution [221 . 

m 


lim <F(z,£) = 

£—>•00 (7 Z 


(70) 


This is derived as the stationary solution of the Fokker- 
Planck equation for the process d X( t ) = <rX(t)dW (t)+dt 
with X(0) = 0. It was shown [22[ that for any t > 0, 
X(t) has the same probability distribution as A(t). For 
intermediate values of t, analytical expressions are also 
available, although they are rather involved, see [Hi HU 
and the references in these papers. The process leading 
to the distribution m is a particular case of a class of 
Markov processes which admit stationary solutions and 
were studied in 24\ and [U HU . 

The probability distribution of y = x(t)/x( 0) 

of the continuous time process has support y : (l,oo). 
In the very large time limit it approaches the stationary 
distribution 


lim ^(y,t) = 


2 r 


a 2 y log 2 y 


-£ <J A log y m 


(71) 


This gives the asymptotic t —> oo probability distribution 
of the random multiplicative process 0 in the continu¬ 
ous time limit under the scaling p = rr. This distribution 
has a power-like tail ~ 0(l/y), up to the slower varying 
logarithmic and exponential factors. 

The expectation value (x(t)) of the solution of the con¬ 
tinuous time diffusion (1551) is infinite under both t —> 0 
and t —> oo limiting distributions. For t —> 0 the variable 
x(t ) is approximative^ log-normal 


x(t) = a:(0) exp(rfe CT \ / *^ A "-3 fT2t ) (72) 


where X is a Gaussian random variable with mean zero 
and variance 1. The expectation of this random variable 
diverges, for arbitrarily small time. A similar divergence 
occurs for the t —> +oo limiting distribution CD- By the 
same argument, all positive integer moments of x(t) are 
infinite. 


VII. SUMMARY AND DISCUSSION 


x{t) 



dse aW{s) -^ 2s ^ . 


( 68 ) 


We considered in this paper the long-run growth rate 
of the average value ( x n ) of a discrete time random multi- 










11 


plicative process m driven by the exponential of a Brow¬ 
nian motion. The study of the n —> oo limit is consider¬ 
ably aided by the equivalence of this average value with 
the partition function of a one-dimensional lattice gas 
with n sites and attractive interaction fl5| . The n —> oo 
limit corresponds to the thermodynamical limit of the 
lattice gas, and can be studied using classical statistical 
mechanics methods. 

The main results of the study can be summarized by 
the following properties of the Lyapunov exponent A: 

i) the Lyapunov exponent has a scaling property in 
the n —> oo limit. In this limit it approaches a function 
A(p, T), with T = 1/P where /3 = i a 1 2 3 4 5 6 t n n is a combina¬ 
tion of the parameters of the process ©, which plays the 
role of the inverse temperature in the equivalent lattice 
gas. The Lyapunov exponent is related to the lattice gas 
pressure p as A = 1 /Tp. 

ii) the functional dependence of A(p, T) for p < pc dis¬ 
plays non-analyticity typical of a first-order phase tran¬ 
sition. The Lyapunov exponent A (p, T) is a continuous 
function of T at fixed p , but its derivative 9j’A(p, T) is 
discontinuous at a transition point T tr (p) provided that 
p < pc- This behavior is related to a phase transition 
in the equivalent lattice gas. The qualitative features of 
this transition are reproduced to a good approximation in 
terms of a mean-field theory and van der Waals equation 
of state [lj. 

The results of this paper can be extended to the Lya¬ 
punov exponents associated with the higher positive mo¬ 
ments of the state variable (( x n ) p } with p > 2. The 
numerical study in jl5[ showed that these moments have 
explosive behavior at certain values of the model param¬ 
eters, similar to that observed for the average value. The 
analogy with the lattice gas can be extended also to these 
moments, which leads to the same scaling property as 
noted for the first moment in the n —> oo limit. 

The model CD can be generalized by replacing the 
standard Brownian motion W* with an arbitrary Gaus¬ 
sian stochastic process sampled on a set of equidistant 


times ti . The equivalence relation with a one-dimensional 
lattice gas © can be extended to such models by re¬ 
placing the two-body interaction energy in Eq. Q with 
£ij = --^cov(Zi 1 Zj). The condition for the existence 
of the Lyapunov exponent is the same as the condition 
for the existence of the thermodynamical limit in the 
equivalent lattice gas. This requires that the covariance 
cov( Z \,, Zj) falls off sufficiently fast with U — j\ such that 
the sum in CD converges for any i [171 ]. However, in 
order for a phase transition to be present, the covariance 
function cov(Zj, Zj) cannot fall off too fast, and the inter¬ 
actions of the equivalent lattice gas must be sufficiently 
long-ranged [29[. Sufficient conditions for the presence 
of a phase transition in a one-dimensional lattice gas are 
given in [29j. 

Assuming translation invariant interactions, the con¬ 
vergence of the sum JT \i — j\£ij has been shown to im¬ 
ply the absence of a phase transition EM. In par¬ 
ticular, this implies that taking X(t) to be a stationary 
Ornstein-Uhlenbeck process for which the covariance de¬ 
creases exponentially co v(X(t)X(s)) = does 

not produce a phase transition unless the mean-reversion 
parameter 7 is scaled in an appropriate way to zero, by 
taking the Kac limit [3lj. Numerical simulations of such 
a model in E3 showed that the average ( x n ) has a very 
rapid explosion with er at a certain point. According to 
the above argument this effect is expected to be smoothed 
out in the n —> 00 limit such that the growth rate of ( x n } 
is an analytical function of the model parameters, and 
no phase transition is present in the analog lattice gas. 
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